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Nonequilibrium states induced by an applied bias voltage (V) and the corresponding current- 
voltage characteristics of one-dimensional models describing band and Mott insulators are investi- 
gated theoretically by using nonequilibrium Green's functions. We attach the models to metallic 
electrodes whose effects are incorporated into the self-energy. Modulation of the electron density and 
the scalar potential coming from the additional long-range interaction are calculated self-consistently 
within the Hartree approximation. For both models of band and Mott insulators with length Lc, 
the bias voltage induces a breakdown of the insulating state, whose threshold shows a crossover de- 
pending on Lc- It is determined basically by the bias Vth ~ A for Lc smaller than the correlation 
length £ = W/ A where W denotes the bandwidth and A the energy gap. For systems with Lc ^> £, 
the threshold is governed by the electric field, Vth/ Lc, which is consistent with a Landau- Zener- type 
breakdown, Vth/ Lc oc A 2 /W . We demonstrate that the spatial dependence of the scalar potential 
is crucially important for this crossover by showing the case without the scalar potential, where the 
breakdown occurs at Vth ~ A regardless of the length Lc- 
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I. INTRODUCTION 

Nonlinear conduction in correlated electron systems 
such as one-dimensional Mott insulators^ and two- 
dimensional charge-ordered materials 3 - — has been of 
great interest in the past few decades. They offer in- 
triguing subjects of nonequilibrium physics in condensed 
matter and possibility for novel functions of electronic 
devices. For example, in a typical quasi-one-dimensional 
Mott insulator, Sr^CuOs,- a dielectric breakdown has 
been observed experimentally by applying a strong elec- 
tric field. A dielectric breakdown has been reported 
also in an organic spin-Peierls insulator, K-TCNQ 2 
[TCNQ=tetracyanoquinodimethane]. For another or- 
ganic compound, (BEDT-TTF) (F 2 TCNQ) [BEDT- 
TTF=bis(ethylenedithio)tetrathiafulvalene], which is a 
quasi-one-dimensional Mott insulator, metal-insulator- 
semiconductor field-effect transistor device structures 
have been reported^ where its field-effect characteristics 
are different from those of band insulators^ - — 

So far, theoretical investigations on nonlinear conduc- 
tion for interacting electron systems that are initially in- 
sulating in their equilibrium states have been done by 
several authors*^"— In general, these studies are classi- 
fied into two approaches depending on whether an ex- 
ternal force is written as an electric fiel d 13 i 14 i 16 ~ — or a 
bias voltage j^ ~ 12 i 15 i 2Q ~ — The former approach is to con- 
sider electron systems without electrodes. The electric 
field is usually applied with open boundary condition^ 
or equivalently with periodic boundary condition by us- 
ing a time-dependent magnetic flux j 13 i 16 In one dimen- 
sion, Oka and Aoki studied the Hubbard model under 
a strong electric field by the time-dependent density- 
matrix-renormalization-group method One of their im- 



portant results is that the dielectric breakdown of Mott 
insulators is interpreted as a many-body counterpart of 
the Landau-Zener (LZ) breakdown ] 23 ! 24 This is known 
to describe the breakdown of band insulators where the 
one-particle picture holds. The threshold is given as 
£?th oc A 2 with A being an energy gap. A mean-field ap- 
proach to electric-field-induced insulator-to-metal tran- 
sitions by using Keldysh Green's functions has been re- 
ported in Ref. [ijj. 

The latter approach is to consider an insulator at- 
tached to electrodes. Interface structures must be ex- 
plicitly taken into account in real electric devices. Along 
this line, Okamoto investigated the current- voltage (I-V) 
characteristics of heterostructures that consist of Mott- 
insulator layers sandwiched by metallic leads by combin- 
ing the dynamical mean-field theory with the Keldysh 
Green's function technique ] 2Q i 21 Ajisaka et al studied 
the I-V characteristics of an electron-phonon system cou- 
pled with two reservoirs by a field-theoretical method. 22 
In the studies listed above, Okamoto discussed the I-V 
characteristics in the framework of the LZ breakdown, 
while Ajisaka et al. proposed a different mechanism with 
the threshold bias voltage Vth ~ A. These results seem 
to be inconsistent with each other. 

When we consider a nanostructure in which some ma- 
terial is attached to left and right metallic electrodes, 
the bias voltage V applied to the material is described as 
V = hl — with fiL and hr being the chemical poten- 
tials of the left and right electrodes, respectively. Here we 
assume that the work-function difference at the interfaces 
is absent for simplicity. In this case, one might expect 
that the current flows when some energy levels of the 
material appear in the region between jul and indi- 
cating the threshold governed by the applied bias voltage 
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V t h ~ A. In fact, for the transport in a field-effect tran- 
sistor with a small channel, such an explanation has been 
used frequently^ However, for a bulk insulator with an 
energy gap A, this picture does not hold and should be 
replaced by the LZ mechanism where the threshold is 
determined by the electric field, E t h oc A 2 . This con- 
sideration poses us a question about which parameter 
determines the mechanism of the breakdown. In partic- 
ular, we address the condition for the realization of the 
LZ breakdown in a structure with electrodes. It is also 
important to know the relation between the approaches 
using the structure with electrodes and those which do 
not include them explicitly. We point out that the size of 
an insulator as well as the potential distribution inside it 
determine the nature of the breakdown. 

In this paper, we study one-dimensional band and 
Mott insulators attached to two electrodes (see Fig. 
[T]), using the nonequilibrium Green's function approach 
that has previously been used to discuss the suppres- 
sion of rectification at metal-Mott-insulator interfaces.— 
This approach more naturally describes nonequilibrium 
steady states than the approach based on the time- 
dependent Schrodinger equation because a current os- 
cillation is inevitable in the latter owing to finite-size 
effects. The present method can be easily applied to 
higher-dimensional systems. Preliminary results for the 
I-V characteristics of two-dimensional charge-ordered 
systems are reported in Ref. [26|. 

We show that the applied bias voltage V induces a 
breakdown of band and Mott insulators at zero tempera- 
ture. For both insulators, the threshold shows a crossover 
as a function of the size of the insulating region Lc- For 
systems with Lc smaller than the correlation length £, 
i.e., the characteristic decay length of the wave function 
in the insulator, the breakdown takes place when the bias 
V exceeds the energy gap A. For Lq ^> £, it is governed 
by the electric field, V/Lq, which is consistent with the 
LZ tunneling mechanism. Whether the charge gap is pro- 
duced by the band structure or by the electron-electron 
interaction is irrelevant to the crossover phenomenon. 
We will focus on the spatial modulation of the wave func- 
tions. For the crossover behavior and the deformation of 
the wave functions, the spatial dependence of the scalar 
potential inside the band or Mott insulator is important. 
This is demonstrated in the Appendix by showing the 
case without the scalar potential, where the bias-induced 
transition occurs at ~ A regardless of the length Lc- 



II. MODEL AND FORMULATION 

We consider a one-dimensional insulator that is at- 
tached to semi-infinite metallic electrodes on the left and 
right sides as shown in Fig. [TJ The insulator is referred to 
as the central part, and described by the Hubbard model 
for a Mott insulator and a tight-binding model with al- 
ternating transfer integrals for a band insulator, both at 
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FIG. 1: (Color online) Schematic view of the model. A one- 
dimensional band or Mott insulator in the central part is con- 
nected by left and right electrodes. Solid (open) symbols rep- 
resent the sites in the central part (electrodes). 
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where c\ a (ci a ) denotes the creation (annihilation) oper- 
ator for an electron with spin a at the ith site, = 
cl a Ci a , and rii = + n^. Lc is the total number of 
sites in the central part. The parameter t denotes the 
transfer integral, St its modulation, and U the on-site in- 
teraction. We use t as a unit of energy. For St = and 
U > 0, the first and second terms in Eq. (pQ) become the 
one-dimensional repulsive Hubbard model, whereas they 
describe a band insulator for St ^ and U = 0. The long- 
range Coulomb interaction term with Vij = V p /\i — j\ 
is introduced because it is responsible for the potential 
modulation near the metal-insulator interfaces (i > 1 
and i < Lc)- The Vij term is treated by the Hartree ap- 
proximation, which is equivalent to the introduction of a 
scalar potential that satisfies the Poisson equation. Here 
((ij)) means the summation over pairs of the ith and jth 
sites with i ^ j in the central part (1 < i,j < Lc)- We 
briefly review our formulation 12 below to treat steady 
states with a finite voltage. 

For the metallic electrodes, we consider noninteracting 
electrons. The effects of the left and right (a = L, R) 
electrodes on the central part are then described by the 
retarded self-energies For simplicity, we take the 
wide-band limit so that the self-energies are independent 
of energy. In the present case, they become diagonal 
matrices 2 ^ 

i i 

(^a)ij = ~2^Xoc)ij = —-jlaSii a Sji a , (2) 

where Su a and 5ji a are the Kronecker delta and %l = 1 
(zr = Lc) denotes the site connected with the left (right) 
electrode. We consider the case of = 1r- Within the 
Hartree-Fock approximation for the on-site term in Eq. 
(PQ), the retarded Green's function for spin a is written 
as 

[G r a (e)- 1 } ij =e6 ij -(H{ lFa ) ij , (3) 
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with 



a=L,R 



(4) 



The off-diagonal elements of (-Shfc)^' come from the first 
term of Eq. (JTJ and the diagonal elements are written as 



(H H Fa)ii = il>i + U(n i9 -1/2), 



(5) 



with a = —a. Here, the scalar potential ipi is defined by 
the Hartree approximation to the long-range Coulomb 
interaction as 



ipi = } Vij({nj) - 1) + ai + b, 



(6) 



with a and b being constants. These constants are so 
determined that ifj satisfies the boundary conditions: 



1 \ for i = 1 
^ = { -£fari = L 0l 



(7) 



for the bias voltage V. When V is positive, the left 
electrode has a higher potential for the electrons and the 
current flows from left to right. 12 Here we assume that 
the work- function difference is absent at the interfaces. 
By diagonalizing the complex symmetric matrix H^ F(J , 
the retarded Green's function is obtained as 



"m(»)"m(j) 

e - E a '' 



(8) 



where E a m = e a rn — (e^ and 7^ are real) is the 

eigenvalue of H^ Fa and u^i) is the corresponding right 
eigenvector. 

The electron density is calculated by decomposing it 
into the "equilibrium" and "nonequilibrium" partsi 2 - as 



(n i<7 )=nZ + Yjnt. 



(9) 



The "equilibrium" part is defined by integrating the local 
density of states as 
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where fc( e ) — #(mc — e) is the Fermi distribution function 
with the chemical potential fic at the midpoint of the 
right and left chemical potentials, fic = (l^R + Ml)/2. 
Since jiL = V/2 and jiR = — V/2, we have \±c = 0. 
For V > 0, where the left chemical potential is higher 
than the right, 5nf a (Sn^) is interpreted as the inflow 
(outflow). 

The "nonequilibrium" part of the density Snf^ is ob- 
tained from the "nonequilibrium" part of the lesser 
Green's function: 



de[8G< a {e) 



(11) 



In order to obtain (5G< a (e), we first decompose the lesser 
self-energy in the wide-band limit 2 ^ as in Eq. ([9]): 

E<(e) = i(T L f L (e)+T R f R (e)) 

= E<^( e ) + ^E<«(e), (12) 



with 



and 



where f a (e) 
equation 



E<^(e) = i(r L + r R )fc(e), (13) 

SE< a (e)=ir a [f a (e)-fc(e)], (14) 
= 0(/jL a — e )- Then, we employ the Keldysh 



5G< a (e)=G r (r (e)5i:< a (e)G^(e), 



(15) 



where G"(e) is the Hermitian conjugate of G r a {e). The 
expressions for and 5nf a , with which the numerical 
calculations are carried out, are obtained by substituting 
Eq. © into Eqs. (JTOj and dTTJ. The results are 
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In the above equations, we recover the electron density 
in the equilibrium state without the electrodes if 7^ = 
0, since the bracket in Eq. (fT6|) is reduced to the step 
function and Snf a = 0. 

The current from the left electrode is expressed by the 
"nonequilibrium" part of the density a o 12 i 3Q 

/°° rJf 
x 2; E^ r ^( e ) r « G "( e )n/^) - fR(t)} 



ILIR 
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where we set e = Ti 
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FIG. 2: (Color online) I-V characteristics of one-dimensional 
band insulators for St — 0.025 and St = 0.05. The other 
parameters are Lc — 200, U — 0, V p — 0.1, and jl = Jr = 
0.1. Results for St — are also shown where we set Lc — 200, 
U = 0, V p = 0, and 7 l = 1r = 0.1. 



III. RESULTS 

In this section, we show the results of I-V characteris- 
tics, charge densities, and the spatial dependence of wave 
functions for band and Mott insulators. For both models, 
a breakdown of the insulating state takes place when the 
bias V becomes sufficiently large. The threshold shows a 
crossover behavior as a function of the size of the central 
part Lc-, which indicates the mechanism of the break- 
down changes according to Lc- The profile of the scalar 
potential ipi has crucial importance on the way of the 
breakdown. This is demonstrated in the Appendix by 
showing that, if ipi is artificially set to zero for all z, the 
crossover phenomenon disappears. 



A. Band insulators 

Figure [2] shows the L-V characteristics for band insu- 
lators with St = 0.025 and 0.05. The other parameters 
are L c = 200, U = 0, V p = 0.1, and j L = = 0.1. For 
comparison, we show the results for the regular trans- 
fer integrals (St = 0) with L c 200, U = 0, V p = 0, 
and 7l = 7^ = 0.1. For St = 0, the current J becomes 
nonzero for V 7^ since the central part is metallic. The 
L-V curve has stepwise structures owing to the finite- 
size effect. For V = and St > 0, the central part is a 
band insulator with the energy gap A = ASt. Because 
of the gap, J is suppressed near V = 0. The I-V curves 
for finite St show more complex structures than that for 
St = 0. Apart from the fine structures, J increases al- 
most linearly for large V, which indicates a breakdown 
of the band insulator. 

For finite voltages applied, the charge distributions in 
resistive and conductive states for St = 0.025 are shown 
in Figs. [3] and HI respectively. In Fig. [3fa), the elec- 
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FIG. 3: (Color online) (a) Electron density, (rii), "equilib- 
rium" part, n® q , and (b) "nonequilibrium" parts Snf and Snf 
for St = 0.025, L c = 200, U = 0, V p = 0.1, lL = Jr = 0.1, 
and V — 0.3. The scalar potential ipi is shown in the inset of 
fa). 



tron density (rii) = (n^) + (n^) and its "equilibrium" 



part n i 



eq 



cq 



for V = 0.3 are shown. A 2fe 



oscillation in the charge distribution is induced by the 
boundaries^ For all z, (n;) is almost unity as in the equi- 
librium case ((fii) = n^ q = 1 for V = 0). The electron 
density (rii) is basically unchanged by the bias voltage 
V when J is small. The scalar potential ^ has a lin- 
ear dependence on i as shown in the inset of Fig. [3ja). 
This is because the long-range interaction term in Eq. 
([6]) is small for (rii) ~ 1? so that ipi is determined only 
by the boundary conditions. The "equilibrium" part n^ q , 
on the other hand, deviates from unity near the left and 
right electrodes, where the deviation is canceled by the 
"nonequilibrium" parts Snf , as shown in Fig. [31(b). The 
quantity Sn\ have nonnegative values for all i because 
electrons come in from the left electrode. Although Sn\ 
is large near the left electrode, it decays as i increases. 
On the other hand, Snf" have nonpositive values for all i 
because electrons go out to the right electrode. Note that 
Snf = —Sn\ c+1 _ i for 7^ = 7^, no work-function differ- 
ences, and at half filling^ The behaviors of Sn\ and Snf 
indicate that electrons and holes hardly penetrate into 
the central part and the resistive state is maintained. 
The current hardly flows through the central part, be- 
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FIG. 4: (Color online) (a) Electron density, (rii), "equilib- 
rium" part, n® q , and (b) "nonequilibrium" parts Snf and Snf 
for St = 0.025, L c = 200, U = 0, V p = 0.1, ~f L = 1r = 0.1, 
and V — 1.0. The scalar potential ipi is shown in the inset of 
(a). 



cause J in Eq. ([T8]) is determined by the difference be- 
tween the density modulation by the left electrode at the 
right boundary, Sn\ R , and that by the right electrode at 
the left boundary, Sri^. Both terms $n\ R and Sri^ are 
vanishingly small, as shown in Fig. [3fb). 

The charge distribution for V = 1.0, where the sys- 
tem is conductive, is qualitatively different from that for 
V = 0.3 as shown in Fig. HI The spatial dependences 
of n| q and (rii) are nearly the same. They increase al- 
most linearly from left to right except in the vicinities of 
the electrodes, where some oscillatory structure appears. 
The distributions of n^ q and (rii) are understood by that 
of ipi shown in the inset of Fig. [4f a). The electron den- 
sity is higher (lower) on the right (left) half where ipi 
is low (high). This behavior is caused by the electrons 
that move through the system in the conductive phase. 
The profile of ipi shows almost a linear dependence on 
i although a small deviation from the linearity near the 
electrodes is visible in contrast to the resistive phase, 
which is because the charge redistribution ((rii) 7^ 1) is 
easier in the conductive phase. As for the "nonequilib- 
rium" parts of the density shown in Fig. 4(b), Sn\ have 
positive values for all i, while Snf are negative for all 
i because the electrons come in from the left electrode 
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FIG. 5: (Color online) I-V characteristics of one-dimensional 
band insulators for several values of Lc with St = 0.025, U = 
0, V p — 0.1, and jl = Jr = 0.1. The solid lines show the 
function J = aVe~ Vth ^ v which fits to the results. 
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FIG. 6: (Color online) Dependence of the threshold bias volt- 
age Vth and the threshold electric field E t h on the size of the 
central part Lc- The other parameters are the same as in 
Fig. [5] The error bars in the fitting are also shown. 



and go out to the right electrode. A finite current flows 
through the central part: 5n^ R and Snf L are finite and 
have the opposite signs. 

Next, we discuss the breakdown mechanism of band 
insulators. In Fig. [5j we show the I-V curves for different 
sizes of the central parts Lc with St = 0.025. To the 
numerical results, the function, 



J = aVe- v */ y , 



(19) 



is well fitted, where a and Vth are parameters. This 
expression originates from the LZ tunneling mechanism 
through which the insulator breaks down with the thresh- 
old voltage Vth*^ For V < Vth, the current J is ex- 
ponentially suppressed due to the energy gap, while it 
increases linearly for V > Vth- When the central part 
is large, the fitting works well as shown in Fig. [5j so 
that the breakdown is consistent with the LZ tunneling 
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picture, although there exist fine structures in the I-V 
characteristics which come from the discreteness of the 
energy spectrum of the central part. As Lc decreases, 
the structure becomes more prominent. For Lc = 50, 
for example, a deviation from the fitting curve due to the 
stepwise structure becomes large, which indicates the LZ 
mechanism is no longer applicable to small-Lc systems. 

Figure [6] shows Vth determined by fitting Eq. (fT9|) to 
the data for each Lc, together with the corresponding 
electric field E t h = Vth/ Lc- For large Lc, Vth is pro- 
portional to Lc, so that E t h becomes a constant. In 
band insulators, the LZ breakdown is known to be in- 
duced by the applied electric field j 23 i 24 Since the one- 
particle picture holds in band insulators, this breakdown 
can be analyzed as a usual interband tunneling problem 
and the threshold electric field becomes E t h oc A 2 /W£^ 
The breakdown occurs when the energy gain by displac- 
ing an electron with charge — e in an electric field E by 
the distance £ = W/A, eE^ overcomes the energy gap 
A. Here, W ~ 4 is the bandwidth and £ ~ 40. We 
have obtained the threshold E t h = 0.0022, which is com- 
parable with the value obtained by the LZ formula ^ 23 ! 24 
(A/2) 2 /v = 0.00125 with v = 2. In short, the threshold 
is governed by the electric field. 

When the central part is small, the fitting to the I- 
V curve becomes worse because the finite-size effect be- 
comes severe. The LZ mechanism is not suitable for un- 
derstanding this breakdown. In this case, another mech- 
anism, in which the threshold is determined by the bias 
voltage, is more appropriate for the following reason. As 
Lc decreases, it eventually becomes smaller than the cor- 
relation length The tunneling occurs when the en- 
ergy gain by displacing an electron by the distance Lc, 
eELc, overcomes the energy gap A. This indicates that 
the mechanism of the breakdown continuously changes 
around Lc ~ £ as a function of Lc- When V exceeds A, 
some energy levels of the central part come in between 
fiL and fiR. For Lc < £, the wave functions do not 
fully decay in the system: the electron injected from the 
left electrode with energy higher than A/2 can reach the 
right electrode through these levels so that the current 
flows. This can be clearly seen in Fig. [5] for Lc = 50 
where the gap is A ~ 0.2 due to the finite-size effect. 
The I-V curve shows an abrupt increase at V ~ A be- 
cause fiL exceeds the lowest unoccupied energy level of 
the central part. Each stepwise increase in the I-V char- 
acteristics corresponds to the increase in the number of 
energy levels located between hl and hr. 

We have numerically confirmed that the results are 
qualitatively unchanged even if the long-range Coulomb 
interaction strength V p and the system-electrode coupling 
strength (= Jr) are varied. Thus, the threshold shows 
a crossover as a function of Lc- When Lc > £, the LZ- 
type breakdown occurs and the threshold is governed by 
the electric field. For Lc < £, on the other hand, the 
current flows when V exceeds the energy gap A. 

It is noted that the spatial dependence of ipi is impor- 
tant for the realization of the field-induced breakdown as 




FIG. 7: (Color online) |^(*i«)^m(^)| 2 witn G =t plotted 
against for V = 0, 0.2, 0.4, and 0.6 in the case of St = 
0.025, L c = 400, U = 0, V p = 0.1, and lL = lR = 0.1. For 
V — 0.2, 0.4, and 0.6, the lines are shifted upward by 0.3, 0.6, 
and 0.9, respectively. 



well as the spatial modulation of the wave functions as 
discussed below. For small V, ^ has a linear dependence 
on i throughout the central part because the electrons are 
localized, (rii) ~ 1, so that the effect of the long-range 
interaction on ipi is small. When the system is conduc- 
tive, the charge redistribution occurs near the interfaces 
where a small deviation from the linearity is seen in 
This charge redistribution weakens the electric field on 
the sites away from the interfaces in the central part. 
For comparison, we show in the Appendix the I-V char- 
acteristics that are obtained by artificially setting ipi = 
for all i. This corresponds to a hypothetical case where 
a sufficiently large charge redistribution occurs near the 
interfaces. There is no electric field in the central part: 
a voltage drop occurs only at the interfaces. In this ex- 
treme case, we obtain the threshold bias voltage Vth ~ A 
regardless of the length Lc, which is in contrast to the 
LZ-type behavior for large Lc in Fig. [6j Such a situation 
never occurs in our calculations with ipi and for realis- 
tic parameters. As we will discuss in the next section 
and also in the Appendix, the effect of the spatial profile 
of ipi on the breakdown mechanism of Mott insulators 
is basically the same as in the case of band insulators. 
Thus, the model without ipi is inappropriate for realistic 
insulators. 

In discussing the breakdown for Lc ^> £, the spatial 
dependences of the wave functions u^i) are crucial as 
explained above. In Fig. [7J we show l^fe^mfe)! 2 
as a function of the real part of the one-particle energy 
ej^ for several values of V in the case of Lc = 400. 
The behavior of this quantity for a =1 is the same. It 
shows whether a given one-particle state contributes to 
the current J. Note that J is obtained by integrating 
\[ G a( e )}iLi R \ 2 over < e < Vl- Since e - e a m ap- 
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pears in the denominator for [G£(e)]^- [Eq. J5J)], the 
one-particle state m with finite \u^ n (iR)u < ^ n (iL)\ 2 in the 
interval — V/2 = fiR < < hl = V/2 gives a large 
contribution to J. This quantity directly shows whether 
the one-particle state m is localized or delocalized be- 
cause it comes from the product of the amplitudes of 
the wave function at the two interfaces %l and %r. If 
\ u m{^^ u rn{^L)\ 2 i s large, the state has finite amplitudes 
at both sides of the central part, so that it is delocalized. 
If l^^iO^mfe)! 2 is small, on the other hand, the state 
has a small amplitude at either of the interfaces. 
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For V = 0, 



{iR)ul n {iL)\ 2 shows two bands that 



0.005 



correspond to the conduction and valence bands in the 
band insulator. Because of the energy gap, no state 
exists in the region -A/2 < < A/2 for V = 
so that the current does not flow at least for V < A. 
When V = 0.2, several states appear in the region 
—A/2 < ej^ < A/2, corresponding to the leakage of one- 
particle states from the electrodes to the central part. 
However, these states do not contribute to the current 
because their l^fe^mfe)! 2 are vanishingly small as 
shown in Fig. [TJ Note that the line for each 7 / 
is shifted upward by |y. For V — 0.4, the number of 

states around = with small \u\ n {iR)u\ n {iL)\ 1 in- FIG. 8: (Color online) Spatial dependence of one-particle 




400 



creases. The energy range where these localized states 
appear becomes wider as V increases. Consequently, the 
delocalized states that contribute to the current depart 
from the region fiR < e < Therefore, the current 
does not flow even if V barely exceeds the gap. As we 
show in the Appendix, the localized states do not appear 
if we set ipi = for all i. It is crucial to take the spatial 
dependence of ijji into account to obtain the modulation 
of the wave functions. 

Figure [8] shows the spatial dependences of the squares 
of the absolute values of the wave functions | 2 for 

V = 0.4 and several m whose are located at the po- 
sitions indicated by the arrows in Fig. [3 Here m is so 
labeled that < e^, for 1 < m < m! < L c = 400. Fig- 
ure [SKa) shows the one-particle state in the lower band 
with = 0.0348 (m = 201) which is inside the gap for 

V = 0. This state is localized on the left half of the cen- 
tral part. The reason is as follows. The scalar potential is 
high (low) near the left (right) electrode. Within each of 
the conduction and valence bands, the state whose weight 
is large near the left (right) electrode has a higher (lower) 
energy than others. In the present case, the indexes m 
for the valence band are m = 206, 203, 201, 199, . . ., 3, 2, 
1, while those for the conduction band are m = 195, 198, 
200, 202, . . ., 398, 399, 400. As m is lowered, the wave 
function of the one-particle state is generally extended to 
a wider region and its largest amplitude is shifted to the 
right, as shown in Fig. [8f b) for the case of ej^ = —0.157 
(m = 191). As m is lowered further, e.g., for ej^ = —0.400 
(m = 175) in Fig. [8fc), the one-particle state is delocal- 
ized to reach the right electrode. Then, its wave function 
has large amplitudes near both electrodes. 

In order to overview the behaviors of the one-particle 
states, we show the contour map of one-particle states 



states f° r ( a ) 171 — 201, (b) rn = 191, and (c) m = 175 

in the case of V = 0.4. The other parameters are the same as 
in Fig. [71 The corresponding are indicated by the arrows 
in Fig. [3 



\ u tn(i)\ 2 on tne {h m ) plane for V = 0.4 in Fig. [9] For 
V = 0.4, the states 190 < m < 200 are localized near 
the left electrode, and the states below are delocalized. 
At the bottom of the lower band, the states are local- 
ized near the right electrode because of the low scalar 
potential near the right electrode. 



B. Mott insulators 

In this section, we consider the case where the cen- 
tral part is described by the Hubbard model. Figure 
[TOl shows the I-V characteristics for St = 0, Lc = 200, 
U = 1.5, V p = 0.1, and j L = Jr = 0.1. When V = 0, 
the system is an antiferromagnetic insulator owing to the 
Hartree-Fock approximation. The energy gap A is then 
0.25. In the previous studies on the I-V characteristics of 
metal-Mott-insulator interfaces ^ 10 ! 11 the current was cal- 
culated by solving the time-dependent Schodinger equa- 
tion. It is argued that the results obtained by the time- 
dependent Hartree-Fock approximation for the electron- 
electron interaction are consistent with those obtained 
by exact many-electron wave functions on small systems. 
For example, the suppression of rectification at metal- 
Mott-insulator interfaces is described by both methods^ 
Although the present time-independent Hartree-Fock ap- 
proximation is worse, we expect the present approach 
captures the essential features of nonequilibrium steady 
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FIG. 9: (Color online) Contour map of one-particle states 
|^m(0| 2 m (h m ) plane for V = 0.4. The other parameters 
are the same as in Fig. 

0.05 
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FIG. 10: (Color online) I-V characteristics of one-dimensional 
Hubbard model for St = 0, L c = 200, U = 1.5, V p = 0.1, and 

1L = 1R = 0.1. 



states under the bias voltage. As we increase V, the cur- 
rent begins to flow at V t h ~ 0.8. The breakdown becomes 
a first-order transition due to the Hartree-Fock approxi- 
mation, which is in contrast to the case of band insulators 
in the previous section. 

The charge distributions in resistive and conductive 
phases at finite V are shown in Figs. [TT]and[T2j respec- 
tively. Their overall features are similar to those in band 
insulators. For all i and V < Vth, the electron density 
(rii) is almost unity, which is basically the same as in the 
equilibrium case ((n*) = n^ q = 1 for V = 0). The scalar 
potential ?pi has a linear dependence on i as shown in the 
inset of Fig. [TTl a). As for the "nonequilibrium" parts, 
Sn\ (Snf) is large near the left (right) electrode and de- 



100 

i 



150 



200 




200 



FIG. 11: (Color online) (a) Electron density, (m), "equilib- 
rium" part, n® q , and (b) "nonequilibrium" parts 5nf and Snf 
for St = 0, L c = 200, U = 1.5, V p = 0.1, ~f L = 1r = 0.1, and 
V — 0.5. The scalar potential is shown in the inset of (a). 



cays as i increases (decreases) [Fig. 11(b)]. Electrons 
and holes do not penetrate into the central part so that 
the current does not flow. 

For V = 1.0 > Vth, the spatial dependences of n^ q 
and (ni) are shown in Fig. [T2l a). They increase almost 
linearly from left to right. This reflects the profile of 
the scalar potential ipi that is higher (lower) on the left 
(right) half. The "nonequilibrium" parts of the densi- 
ties, 5n\ and Snf are extended over the whole system 
with small spatial dependence [Fig. 12(b)] ^ which is in 
contrast to the resistive phase. Since Sri]" and Snf are 
finite with opposite signs, a finite current flows through 
the central part. 

In Fig. [13j we show the threshold bias voltage V^h and 
the corresponding electric field E t h = Vth/Lc for the 
first-order transition as a function of Lc- When the cen- 
tral part is small, i.e., Lc ^ 50, Vth is almost a constant 
near the energy gap A c± 0.25. For small Lc, the elec- 
tron injected from the left electrode with energy higher 
than A/2 can reach the right electrode since the correla- 
tion length £ = W/A c± 16 is comparable to Lc- Thus, 
the threshold for Lc < £ is determined by the bias volt- 
age. For Lc ^> £, on the other hand, Vth is proportional 
to Lc, so that the threshold is governed by the electric 
field. In recent theoretical studies, 14 Oka and Aoki have 



9 




200 



200 



FIG. 12: (Color online) (a) Electron density, (m), "equilib- 
rium" part, n® q , and (b) "nonequilibrium" parts 5nf and Snf 
for St = 0, L c = 200, U = 1.5, V p = 0.1, ~f L = Jr = 0.1, and 
V = 1.0. The scalar potential ipi is shown in the inset of (a). 



Vth 0.8 




100 200 300 400 

L C 

FIG. 13: (Color online) Dependence of the threshold bias 
voltage Vth and the threshold electric field E t h on the size of 
the central part Lc- The other parameters are the same as 
in Fig. [10] 



proposed that the LZ breakdown occurs also in Mott in- 
sulators by applying the time-dependent density-matrix- 
renormalization-group method to the one-dimensional 
Hubbard model under an electric field with open bound- 
ary condition. In our calculations, the scalar potential 
ipi is linearly increasing with i in the resistive phase as 



shown in Fig. [TTl(a). which means that the electrons feel 
a uniform electric field in the central part. Therefore, 
our model describes the LZ breakdown as in the open 
Hubbard chain as long as the electrodes do not affect 
the nature of the breakdown for Lc ^> £• In fact, the 
threshold £? t h is about 0.0033, which is comparable to 
the LZ value, 14 (A/2) 2 /v = 0.0078, with v = 2. Thus, 
the threshold shows a crossover as a function of Lc as in 
the case of band insulators. 

The LZ breakdown is explained as before by compar- 
ing the charge gap A and the work which is done by the 
electric field on an electron moving over the correlation 
length £. If the work eE^ exceeds A, the electron in 
the lower band may go over to the upper band so that 
the current flows. According to the results by Oka and 
Aokij 1 ^ this consideration is applicable to Mott insula- 
tors where the correlation effects are important. Thus, 
we expect that the results obtained by the Hartree-Fock 
approximation are qualitatively unchanged even if we 
take account of the electron correlation. It is well known 
that the Hartree-Fock theory overestimates the charge 
gap A. It predicts the antiferromagnetic spin ordering 
which is actually destroyed if quantum fluctuations are 
appropriately taken into account. However, the overes- 
timated A will alter the threshold only quantitatively. 
We also expect that the spin ordering does not essen- 
tially affect the breakdown itself since only the charge de- 
grees of freedom are relevant to the mechanism. We note 
that the time-dependent density-matrix-renormalization- 
group method has been applied to a Mott insulator with 
electrodes very recently^ and that the L-V characteris- 
tics have been consistently explained by the LZ tunneling 
mechanism. However, it is also shown that physical quan- 
tities such as the spin structure factor and the double oc- 
cupancy do not reach a stationary state in the accessible 
time window. Therefore, a direct description of nonequi- 
librium steady states as in the present study is important 
to deepen our understanding of the breakdown. 

Since the crossover behavior is obtained for both band 
and Mott insulators, the electron-electron interaction is 
not responsible for the phenomenon. We emphasize that 
the spatial profile of ipi is important for the realization 
of the LZ breakdown. In the Appendix, this is demon- 
strated for the Mott insulator by showing the I-V char- 
acteristics that are obtained by artificially setting ipi = 
for all i. In this case, no electric field exists in the cen- 
tral part so that the modification of the wave functions 
does not occur. The I-V curves do not show any Lc 
dependence apart from the fine structures coming from 
the discreteness of the energy spectrum. We obtain the 
threshold bias voltage V t h ~ A regardless of the length 
Lc as in band insulators. 

In Fig. [TH we show {u^i^u^iL)] 2 as a function 
of the real part of the one-particle energy for sev- 
eral values of V, where each line for V ^ is shifted 
upward by V. The behavior of this quantity for a =1 
is the same. For V = 0, no state exists in the re- 
gion —A/2 < < A/2 since the energy gap opens. 
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FIG. 14: (Color online) |iC(*flXi(*iO| 2 with a =t plotted 
against for F = 0, 0.3, 0.6, and 0.9 in the case of St = 0, 
L c = 200, U = 1.5, V p = 0.1, and lL = lR = 0.1. For 
V — 0.3, 0.6, and 0.9, the lines are shifted upward by 0.3, 0.6, 
and 0.9, respectively. 



3 0.02 




\ u m(^R) u m(^L)\ 2 shows two bands that correspond to 
the upper and lower Hubbard bands. When V = 0.3, 
one-particle states leak from the electrodes to the cen- 
tral part so that several states appear in the region 
— A/2<eJ n <A/2. As in band insulators, these states 
do not have any contributions to the current because 
their |ixJ n (i^)^J n (i^)| 2 are vanishingly small. As we in- 
crease V further, e.g., V = 0.6, the states with van- 
ishingly small l^^iO^mfe)! 2 appear in a wider range 
around 6^ = 0. The number of these localized states 
also increases. The appearance of the localized states 
keeps the central part resistive until V/Lc reaches the 
threshold electric field E t h even if V > A holds. When 
the system is conductive (V > Vth)> the upper and lower 
Hubbard bands are merged into a single metallic band. 
In this case, all the states around ~ are delocalized 
and contribute to the current. 

Figure [15] shows 1^(0 1 2 as a function of i for V = 
0.6. Here m is chosen at 100, 90, and 85, whose ej^ 
are located at the positions indicated by the arrows in 
Fig. [TH In Fig. [T5l a), we show the wave function with 
e m = —0.00395 (m = 100), which is localized on the left 
half of the central part. This state belongs to the lower 
Hubbard band. Since the scalar potential is high (low) 
near the left (right) electrode, the state is located near 
the top of the lower Hubbard band. The one-particle 
states in the band gradually lose their localized nature 
as m is lowered. This can be seen in Fig. fTST b) for 
the case of ej^ = —0.361 (m = 90), where its largest 
amplitude is shifted to the right compared to that of m = 
100. Figure [T5tc) shows the one-particle state for ej^ = 
—0.497 (m = 85), which is completely delocalized. Its 
wave function has large amplitudes near both electrodes. 
The spatial dependences of the one-particle states in the 



FIG. 15: (Color online) Spatial dependence of one-particle 
states |^Jn(0| 2 f° r ( a ) 171 — 100, (b) m = 90, and (c) rn = 85 
in the case of V = 0.6. The other parameters are the same as 
in Fig. 1141 The corresponding are indicated by the arrows 
in Fig. [H 



resistive phase are similar to those in band insulators. 



IV. SUMMARY 

We have investigated the I-V characteristics of the 
one-dimensional band and Mott insulators attached to 
electrodes. A tight binding model with alternating trans- 
fer integrals for the band insulator and the Hubbard 
model for the Mott insulator are studied by using the 
nonequilibrium Green's function method. The applied 
bias voltage induces a breakdown of the insulating state 
to convert into a conductive state for both models. The 
threshold shows a crossover as a function of the size Lc 
of the insulators. For Lc ^ £ = W/A, the breakdown 
occurs at Vth ^ A so that the threshold is governed by 
the bias voltage. For Lc ^> £, the electric field deter- 
mines the threshold, Vth/Lc oc A 2 /W, which is con- 
sistent with the LZ breakdown reported previously ] 13 i 14 
Since the crossover is obtained for both band and Mott 
insulators, the electron-electron interaction is not respon- 
sible for the phenomenon. The profile of the scalar poten- 
tial ipi , which is linearly increasing with i in the resistive 
phase so that the electrons in the central part feel an 
almost uniform electric field, is important for the realiza- 
tion of the LZ breakdown and the crossover behavior. 
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Appendix 

In this appendix, we show the results when we artifi- 
cially set ipi = for all i. The I-V curves of the band in- 
sulator with St = 0.05, U = 0, 7l = 7i? = 0.1, and ipi = 
for L c = 100, 200, and 400 are shown in Fig. M It is 
apparent that the breakdown occurs at Vth ~ A regard- 
less of the length Lc, which is consistent with Ajisaka et 
a\2& This is in contrast to the results in Fig. [6] where Vth 
is proportional to Lc for Lc ^> £ • 

When we fix ipi = 0, one-particle states do not leak 



from the electrodes to the central part since the electric 
field is absent from the central part. In Fig. [T7J we plot 
\ u rn(^R) u m(^L)\ 2 as a function of ej^ for several values 
of V. The results indicate that the one-particle energies 
and the wave functions are not affected by V. This comes 
from the fact that the Hamiltonian in Eq. (1) does not 
depend on (m) for band insulators if we set ipi = for 
all i. Therefore, no localized state appears inside the gap 
for V = 0. In this case, the current begins to flow when 
V merely exceeds Vth ~ A since one-particle states with 
finite l^^iO^mfe)! 2 appear in the region jar < < 
Ml- 

In Fig. [THl we show the I-V curves of the Mott in- 
sulator with St = 0, U = 1.5, 7l = 7i? = 0.1, and 
ipi = for L c = 100, 200, and 400. As in the case 
of band- insulators, the breakdown occurs at Vth ~ A for 
all Lc- Although the breakdown seems to be continuous 
for Lc = 100 and 200, a small discontinuity is evident for 
Lc = 400, which indicates a first-order transition. The 
discontinuity is more obvious for large U as shown in the 
inset of Fig. [[8] for U = 2 with the gap A = 0.68. 
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FIG. 17: (Color online) |iXm(*fl)^m(^)| 2 with <* =t plotted 
against for V = 0, 0.2, and 0.4 in the case of Lc = 200. 
The other parameters are the same as in Fig. 1161 
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FIG. 18: (Color online) I-V characteristics in the case of St = 
0, U = 1.5, 7l = 7# = 0.1, and tpi = for several values of 
Lc- The results for U = 2 and tpi = are also shown in the 
inset. The arrow indicates the location of the gap A. 



